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Abstract 

The advection of a passive scalar by a quenched (frozen) incompressible veloc- 
ity field is studied by extensive high precision numerical simulation and various 
approximation schemes. We show that second order self consistent perturbation 
theory in the absence of helicity perfectly predicts the effective diffusivity of a 
tracer particle in such a field. In the presence of helicity in the flow simulations 
reveal an unexpectedly strong enhancement of the effective diffusivity which is 
highly nonperturbative and is most visible when the bare molecular diffusivity of 
the particle is small. We develop and analyse a series of approximation schemes 
which indicate that this enhancement of the diffusivity is due to a novel second 
order effect whereby the helical component of the field, which does not directly 
renormalize the effective diffusivity, enhances the strength of the non helical part 
of the flow, which in turn renormalizes the molecular diffusivity. We show that this 
renormalization is most important at low bare molecular diffusivity, in agreement 
with the numerical simulations. 



DAMTP-2000-139 
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1 Introduction 



The advection of passive fields subject to molecular diffusion and convection by turbu- 
lent fluid has been extensively studied by both theoretical and computational techniques 
HI, ||, ||, f|, ||, § . By comparing the results of simulation with the theoretical prediction 
for various long-range quantities, the efficacy of the theoretical methods can be tested 
albeit, in somewhat artificial models. The applications to the physics of complex sys- 
tems and engineering are many fold. In practical problems we need to calculate the bulk 
properties of random media given statistical models for the disorder present. In general 
the complexity of these real world problems means that one must resort to approxima- 
tion schemes to calculate these large scale bulk properties. It is therefore essential to 
verify various methods of analysis on model problems before one can be confident that 
these or similar methods can be applied to more realistic systems. The success of an 
approach depends on whether the approximation preserves the essence of the physical 
mechanism responsible for determining the long-range parameters of the advection in 
terms of the parameters describing the local characteristics of the flow. In this paper 
we consider advection in a helical Gaussian turbulent flow which was originally studied 
in |Q] . The surprising result observed on the basis of simulation is that the long-range 
effective diffusivity, K e , is greatly enhanced by the presence of the helicity by more than 
a factor of two, the effect being strongest for small molecular diffusivity, kq. In the ab- 
sence of helicity the calculation of K e to two loops in self-consistent perturbation theory 
agrees accurately with the simulation for all kq. However, such an approach predicts 
that even maximal helicity will have only a small effect of the order of 10%. This is 
in stark contrast to the results of simulation. The puzzle is to explain these results for 
what is a relatively simply posed system. A successful theoretical approach will involve 
infinite resummations of contributions and it is in this sense that the enhancement is 
non-perturbative. 

In this paper we discuss a possible resolution of the conflict between theory and 
simulation by using various methods to identify the low-wavenumber effective theory 
governing the diffusive dispersal of particles advected in the turbulent flow when helicity 
is present. The derivation of the effective theory is guided by the renormalization group 
(RG) idea that the Green function at low wave-number is, in some approximation, 
the solution to an effective second-order differential equation whose parameters are 
determined self-consistently in terms of the original or 'bare' defining the model. The 
effect of helicity in the flow causes the turbulent velocity field, u(x,t), to be additively 
renormalized by a term proportional to the vorticity, ld = V x u. The coefficient of 
proportionality is a is a pseudo-scalar which is generated by the axial-vector nature of 
the helical flow and so depends on the helicity h, defined in terms of u by 



where (•) denotes the ensemble average over the random velocity field. In our model 
the magnitude of h is measured by a parameter A, < A < 1, and the results are given 
in terms of A. The usual perturbative result for the dependence of n e on A is that n e 
is a series in A 2 for all values of kq. This is self-evident since the magnitude of K e is 
independent of the sign of A. The simulation is seemingly consistent with this fact for 
A < 0.2 at kq = but is not well fitted by any simple approach, and for larger A the 




(1.1) 
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curve lies far higher than the naive calculation. We discuss an improved self-consistent 
scheme which expresses the Green function and vertex functions as solutions to integral 
equations which are solved in a low-wavenumber approximation. This method leads 
to a strong enhancement of K e for increasing A and, as such, is a good indication that 
we are on the right track. However, for small A the effect is paradoxically too strong, 
leading to a non-analytic dependence of n e on A which is predicted to be n e ~ A 2//3 in the 
one-loop case. This is possibly due to the approximation made in obtaining the solution 
but it is a complex matter to ascertain whether this is so. An alternative approach is 
to use the functional Hartree-Fock method which leads to an integral equation for the 
Green function self-energy as a function of wave-number. The result of this method for 
K e (A) is better behaved at small A but the predicted enhancement is not big enough 
and does not fit the simulation data. In general, the effect is most pronounced for small 
Ko and empirically from our simulation we find that the results distinguish the regions 
kq <C 0.2 and kq 3> 0.2. There is a pronounced dip in K e vs Ko at ~ 0.2 for A = 1. 
This dip is not predicted by either of the methods mentioned so far. 

We also present a renormalization group approach which shows a mechanism for the 
enhancing effect of helicity on K e . The renormalization group is normally most useful 
for computing anomalous exponents since they are generally independent of much of the 
details defining the model: the idea of universality. It is much more difficult to control 
a standard RG analysis if it is used to calculate the coefficients of scaling behaviour, 
i.e., observables like K e . However, in ref. (Tj we reported on a successful use of the RG 
in predicting K e for gradient flows and believe that an RG analysis can generally give 
a strong indication of the kind of mechanism which influences the size of parameters 
controlling the large-scale characteristics of advection. In this paper we show that the 
flow at large wave vector can strongly enhance n e when kq is small. In particular, this 
approach does provide a mechanism for the dip observed in K e vs k,q at kq ~ 0.2 for 
A = 1. 

In section || the model and the formalism are reviewed; in section || the perturbation 
theory is briefly described; in section |] the self-consistent integral equations for the 
Green function and vertex functions are derived to one-loop and the small wave vector 
approximation for K e is derived; in section |5| the functional Hartree Fock method is 
examined; in section ^ the renormalization group approach is explained and in section 
the conclusions are presented. 

2 The Model and Formalism 

In [||] the problem of a passive scalar advected by an incompressible turbulent flow 
with a molecular diffusivity was studied. The turbulent fluid velocity field, u(x, t), was 
described by its statistical properties which were assumed to be Gaussian and so fully 
determined by the velocity auto-correlation function. In the original study the flow 
was time-dependent, but since the enhancement of n e by helicity in the flow is present 
also for time-independent flows we assume here, for simplicity, a time-independent flow 
(i.e. quenched or frozen turbulence) for which the auto-correlation function can be 
expressed in the following form: 

^(xK(x')> = J ^ e *-(x-x0^. (k) . {21) 
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The ensemble of velocity fields was taken to be homogeneous and isotropic and so 
for incompressible fluids .Fy(k) can be written as 

Fij(k) = $(k)(fc 2 <% - hkj) + y(k)ie imj k m , (2.2) 

where ^ represents the presence of helicity in the flow. In [|j it was assumed that <3? 
and \E' took the factorized forms: 

*(k) = &^-A*E(k) 

tf(k) = ^^-A 2 kE{k)sm 2il) , (2.3) 
3 

where A is chosen so that 

' dkE(k) = 1 , (u • u) = «g , (2.4) 



and where uq is the r.m.s. velocity. Choosing the angle tjj to be ^-independent means 
that the helicity is of equal strength at all wave vectors. The helicity parameter, h, has 
been defined in eqn. (1.1) and with the definitions in eqn. (2.3), we find 

h = ^A 2 (k 3 )sm2ip , (2.5) 

where (k 3 ) is the expectation value with respect to the distribution E{k). The passive 
scalar field B(x, t) is advected according to the equation 

^ = K0 V 2 G- V-(u0) , (2.6) 
at 

and the effective, or long-range, diffusivity, K e , is defined by 

(x • x)(t) = (J d 3 x x • x 6(x, t)) , 

= 6M + O(t ) as t -> cc , (2.7) 

where is normalized to unity: 

r d 3 x 0(x,t) = 1 . (2.8) 



For the purposes of numerical simulation a particular member of the velocity-field 
ensemble is then realized by ^, |j 



N 



u(x) = A { U n cos ^ ~ X n A kn sin V^) A k n cos (k n • x) 

n=l 

k n sin A k n sin 

(k n • x) } , (2.9) 

where the vectors £ n and Xn are distributed uniformly and independently over the unit 
sphere and the wave vector k n is distributed according to the distribution E(k). For 
N sufficiently large the central limit theorem guarantees that u(x) is Gaussian up to 
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0(1/ N) corrections. We have used N = 64 for which these effects are sufficiently small 
for our purposes. 

To simulate the evolution of the scalar field 0(x, t) we integrate numerically the 
stochastic equation for the evolution of a particle with path x(t) given by 

x(i) = u(x(i)) + rj(t) , (2.10) 

where r/(t) is a Gaussian random variable with (rj(t)) = and (r](t)-ri(t')) = 2Ko5(t—t'). 
The resulting probability distribution for particle position x(i) is then 0(x, t) with the 
initial condition 0(x, 0) = 5(x). 

The discrete form of eqn. fl2.10|) suitable for numerical integration is: 

x n+1 -x n = u(x n ) At + (2rc Ai)2 e n , (2.11) 

where e n is a Gaussian random three-vector of zero mean and unit variance for each 
component. This equation models eqn. (|2TT0| ) correctly to O(At) but the details of 



a third-order Runge-Kutta scheme correct to 0(At 3 ) are given in [Q]. We use this 
third-order scheme in our numerical simulation. 

The effective diffusivity, K e , is then computed from the ensemble of paths by 

1 M 

" x(t)) paths = Jim^-E^W-x^W, 

a=l 

= 6n e t + 0(f) as t -> 00 . (2.12) 

Here M is the total number of paths averaged over and (a) label the member of the 
ensemble of paths. In practice M is finite but large enough to give an estimate for 
K e with small error. In addition t must be large enough so that the path evolution 
is in the asymptotic regime where the evolution can be suitably described in terms of 
long range effective, or "renormalized" quantities. That t is large enough is tested by 
ensuring that the estimate for K e is independent of t within statistical errors. 

3 Perturbation Theory 

The perturbative approach to solving eqn. ([O]) is well known || |!| |(| and we only 
summarize here the necessary results. 

Since we are interested in the effective parameters governing the evolution of the 
distribution 0(x, t), we study the related Green function G(x) which satisfies 

k V 2 G -u-VG = - <5(x) , (3.1) 

where the incompressibility of u has been used. A perturbation series in u/fco for G(k) 
can be generated by iterating the formal solution to eqn. ( p.l| ) in Fourier space: 

The Green function averaged over the velocity ensemble, (G(k)), can then be written 

as 

(6(k » - ^hw) ■ (3 - 3) 
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where the averaging over the velocity ensemble is done using Wicks theorem to give a 
diagrammatic expansion and E(k) is given by one particle-irreducible diagrams. The 
asymptotic behaviour in eqn. ( |2.12| ) implies that the small k behaviour of (G) is given 
by 

The Feynman rules for the diagrammatic perturbation expansion are: 

(i) Wavevector is conserved at each vertex; 

(ii) Each full line carries a factor of 1/k 2 ; 

(iii) Wavevector is integrated around closed loops with a factor dq/ (27r) 3 ; 

(iv) The primitive vertex Vj(k',k), whose diagrammatic representation is shown in 
figure ||a, is given by Vj(k', k) = i k'i ; 

(v) Each internal dotted line carries a factor 

Fy(k) = $(k)(fc 2 ^- - kikj) + ^{k)ie imj k m . (3.5) 

In what follows we use the explicit spectra 

•00 = ^f^- k> °. 

tf(k) = Afc*(k), (3.6) 

where A = sin 2ip . 

The simple two-loop calculation for n e gives the result 

flu 2 u 4 \ 

« e = «o 1 + Q^fa + (0.0059A 2 - 0.00884)^5 . (3.7) 

The diagrams contributing to this order are shown in figure ^. The two-loop integrals 
were done numerically. 

The effect of helicity is not seen until second order. This is evident from the ex- 
plicit expressions for the diagrams but is also easily understood because the effect of 
helicity on n e cannot depend on the sign of A. Hence, the graphs with a non-vanishing 
contribution from helicity must contain an even number of internal velocity correlators 
(dotted lines). 

Clearly, this approach is not applicable to the limit kq — > in which we are interested 
but a self-consistent method will allow the model to be perturbatively analyzed in this 
limit and this is described in the next section. 



4 Self-Consistent Methods 

A self-consistent treatment performs a resummation of an infinite subset of diagrams 
which gives a continuation of the perturbation theory beyond its normal radius of 
convergence. The approach is not unique but depends on how the effective low energy 
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theory is parametrized and which quantities are treated self-consistently. A successful 
result will depend on how well the method captures the dominant effects in this way. 

We first discuss the simplest approach which treats only K e self-consistently. At 
two-loops this gives an excellent fit for n e when A = but fails for A / 0. We then 
generalize the method and show that we can qualitatively explain the large enhancement 
in K e due to helicity although the approach is still quantitatively deficient. Further 
generalizations are discussed but have not yet been carried out. 



4.1 Self-Consistency in n e 

To generate the self-consistent perturbation series in K e the eqn. for G(x) is 

formally rearranged to become 

k £ V 2 G-AkV 2 G-u-VG = -<5(x), (4.1.1) 

where An = n e — kq . The second term is a counter-term which is included as part 
of the perturbation. It is formally of first order in the expansion parameter which 
allows the expansion for K e to be constructed to a consistent order. The self-consistent 
perturbation series is generated by iterating 

^ (k) = i + ^? / Wf ( <(k " q) ' Hq) " Ak 6{q) (k " q)2 ) 6(k " q) ' (4 ' L2) 

This equation can be rewritten as an equation for Ak and its diagrammatic represen- 
tation up to two loops is shown in figure |3](a). Since K e is not renormalized from the 
tree-level value we have the self-consistency condition 

A S(k) | k=o = o . (4.1.3) 
To A-th order in u^/ K^k^ it is always possible to write this condition in the form 

K e = K + K e F N (K e ,\) , 
N / 2 \ n 

F N ( Ke ,X) = £<?«(A) o . (4L4) 

n=l \ K e K 0/ 

^From now on we set uq = 1, ko = 1 and, using the velocity-field spectrum given 
in section g we find the two-loop self-consistent expression for n e 

( Ke - K )(l - = Kef^ + \ fo.0059A 2 - 0.00884]) . (4.1.5) 



This result can be re-expressed in the form of eqn. ( 4.1.4 ) to become 

/ 1 0.0035 + 0.0059A 2 \ , J1fl . 
Ke = n + Ke ^ + - t j . (4.1.6) 

We show the two- loop self-consistent prediction for n e compared with data in figures [7] - 
ll. In figures 0, ^ and |^ K e is plotted against kq for fixed I = 0.0, 0.4, l.O and in figures 



10 and |ll] K e is plotted against A for fixed kq = 0.0, 0.2. As should be expected, we see 
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from figures |8] and |9] the agreement between theory and simulation is acceptable for kq 
large enough. This is simply because the large molecular diffusivity swamps all other 
effects. However, there is a large disagreement for small kq which is most marked for 
kq = 0. The prediction for n e behaves like 0(X 2 ) and for n$ changes from n e = 0.3697 
at A = to K e = 0.4090 at A = 1: an increase of 10%. In contrast, the simulation gives 
K e = 0.3705(1) and K e = 0.8018(7) respectively at these two values of A: an increase of 
more than a factor of two. From the simulation for kq small enough we find that n e 
as a function of A is strongly in disagreement with the slow polynomial behaviour in 
A predicted by self-consistent perturbation theory. This effect was first observed in j|] 
and has remained unexplained. 

In addition, in figure [9| we observe a marked dip in the data at fixed A = 1 for 
K e versus kq at about kq = 0.2. The major feature is that K e rises rapidly with A at 
kq = whereas the effect for kq > 0.2 is much less strong: the dip is not a lowering 
of the curve as A increases at kq = 0.2 but rather a rapid rise with A at kq = 0. The 
self-consistent prediction of this section does not predict a dip of any kind. 



4.2 A more general approach 

In this section we propose an explanation of the enhancement of K e by helicity in the 
flow. The technique is presented in detail at the one-loop level and the extension 
two-loop is then given. 

The philosophy is to suggest an effective, low wavevector diffusion equation obeyed 
by the smoothed distribution function. Because the wave vector is small it is assumed 
that the equation can be limited to at most two spatial derivatives. The shortcomings 
of this assertion are discussed later. We propose the equation 

^ = Ko V 2 e-a R u-VQ- p R u-VQ , (4.2.1) 

where u is the vorticity and ocr and (3r are coupling constants which must be de- 
termined self-consistently. The bare values of these couplings, ao and 0o, define the 
original diffusion equation. From eqn. (^) we see that in our case ao = 1 and /J = 
but the analysis may be applied for general values of ao and (3$. Eq. ( 4.2.1]) is equiva- 



lent to renormalizing the velocity field to ur = a^u + 0rU). It will turn out that a is 
not renormalized and so the self-consistency conditions are applied only to determine 
K e and (3r. To this end the equation for the G(x) is taken to be 

K e V 2 G — An V 2 G — ao u • VG — [3ru • VG + A/3 u • VG = - <J(x) , (4.2.2) 

where, as before, A and A/3 = (3r — (3q . The rules for perturbation theory 

are the same as given in section || with the additional rules 

(vi) The primitive vertex Wj (k , k) associated with the vorticity and whose diagram- 
matic representation is shown in figure |l]b, is given by Wj(k',k) = (k' A k)j; 

(vii) For each vertex of type Vj a factor of a# and for each vertex of type Wj a factor 

Of 0R. 
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The self-consistent equations are given by setting the next renormalizations of n e and 
(3r to zero in perturbation theory. This gives two equations which simultaneously 
determine K e , or and (3r in terms of the bare parameters kq , a$ an d Po ■ It is convenient 
to define the general vertex Uj(k',k) of the form 

U;(k', k) = iV(k', k) k'i + W(k', k) (k' A k)j , (4.2.3) 

where the form- factors V and W are scalar functions of k and k'. The bare vertex 
U° is defined by V° = a , W° = fa. There is no independent form-factor coefficient 
proportional to kj in this expansion since the velocity field is incompressible. The 
diagrammatic representation of U° is shown in figures |l]c and |l]d, where the bare vertex 
is represented by an open circle while the renormalized vertex carries additionally an 
inset letter 'R'. Likewise, the general expression for (G(k)) can be defined as 

<G(k)> = 1 (4.2.4) 

k 2 (kq + n(k 2 



To two loops the self-consistent relationships that hold between diagrams are shown in 
figure ||. 

We are interested in the low-wavenumber properties of the theory and we use the 
approximations 

n(k 2 ) = ( Ke - «o) + o(k 2 ) , 

V R = a R + 0(k 2 ,k ,2 ,k-k') , 

W R = R + 0(k 2 ,k /2 ,k-k') . (4.2.5) 

These approximations are consistent with the form considered for the effective equation 
governing the dispersal, eqn. ( [4.2. 1| ), since the renormalized couplings, aR and (3r, are 
given respectively by the coefficients of k j and (k' A k)j in the expansion of Uj. To 
include higher powers of k and k' would, for consistency, require terms involving u and 
u in eqn. ( [4.2. 1| ) with higher powers of V. In principle, a functional self-consistent 
formalism could then be set up for O(k) and Uj(k',k) treating them respectively as a 
vector and matrices in wave-number space. However, this is a further generalization 
that we have not yet pursued because the amount of computing effort required to 
determine them numerically is prohibitive. In the present case we are interested only in 
the low-wavevector behaviour of these functions, and then the self-consistent equations 
are determined by asserting the self-consistency in the limit k',k — > and using eqns. 



(HID 



The calculation is now straightforward but tedious. All loop integrals are approxi- 
mated by their lowest non-zero power in (k, k') and their contribution to the renormal- 
ization of the relevant coupling constant is read off. It suffices to give some examples 
which indicate how the full result is obtained and to this end we first analyze in detail 
the one loop approximation to the self-consistent equations. The first observation is 
that a is not renormalized, i.e., aR = olq. We give one example indicating how this 
comes about. The graphs is labelled by the couplings corresponding to the types of 
vertex it contains and the label is ordered in the same order that they occur in the 
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graph. We consider the contribution shown in figure || to the vertex renormalization 
and will concentrate on the part proportional to k'. The value of this graph is 

rp 2 R f dq ei mn k m q n (k' - q) p k'iFi p (q) 

= "rPrJ (2vr)3 Ke (k-q)^ e (k'-qF " ^ 

The approximations of eqn. ( |4.2.5| ) have been implemented. Only the helical part of 
Fij(q) contributes and we find the result 

a R /3 R X , f dq ei mn k m q n (k' - q) p q *(q) 
Paa k\ 1 V (2vr)3 (k-q)2(k'-q)2 

«|feAy f dq (k-k / q 2 -k-qk / -q)g^(q) 

K 2 lKl J (2vr)3 (k-q) 2 (k'-q) 2 ' 1 

Clearly the contribution to V is O(k-k') and so a is not renormalized. All contributions 
to V are similarly of higher order and the result is that a R = ao . 

The coupling (3 is renormalized when A ^ 0. The calculation follows a similar path 
to that used in the analysis of the renormalization of a. Again we show one calculation 
explicitly and consider the contribution to W(k, k!) by calculating the coefficient of 
(k' A k) in T Qaa : 

t - - a R f dq (k-q);(k / -q)ni ; )n(q) s x 

aaa ~ \lJ (2Ttf (k-q) 2 (k'-q) 2 • 

The contribution proportional to (k' Ak) comes only from the helical part of i^ n (q) and 
so the relevant term is 

a%X f dq e; n pk;k' w q p q m g$(q) 
aaa ~ k\ J (2tt)3 (k-q)2(k'-q) 2 ' 1 ' ^ 

Hence we find the contribution 8(3 R to the renormalization of (3 from T aaa to be 

The renormalization of /9jj is expressed in terms of three integrals: 

J n = ^3 / dqq n <S>(q) n = 1,2,3 . (4.2.11) 
After evaluating all the relevant graphs the self-consistent equations are 
olr - a = , 

Po-Pr + (4A/i + 2a 2 R p R I 2 + a R (3 2 R \h) = 0. (4.2.12) 

The approximate equation for G(k) is given by the equation for £(k) in terms 
of the one-particle irreducible graphs in figure || at one-loop order. Because we are 
using the low-wavenumber approximation this reduces to substituting the expression 
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for the renormalized vertex Uj(k',k) given in eqns. ( 4.2.3j ) and ( [4,2. 12j ) into the one- 
loop diagram for S(k) in figure |5| We analyze the one- loop self-energy graph and keep 
only the term proportional to k 2 . In obvious notation this gives the results 



2a 



R T 



2a R /3 R \ 



(4.2.13) 



Using the spectra in eqn. (|3.6| ) eqns. ( 4.2. 12| ) and (4.2.13) for the one-loop self- 
consistent conditions become: 



Aa = a R — ao 







A/? + ^ = , Ak + ^ 







(4.2.14) 



where 



-3- ( 2a 2 R p R + 2 J - a R 0R A + J- a|A 

lo \ V 7T V 7T , 

~l (a% + A^a R p R X + 3f3 R j . 



(4.2.15) 



^From these equations it is clear that no renormalization occurs if there is no pseudo- 
scalar or axial-vector quantity in the problem: if (5q = A = then the problem reduces 



to the one-loop self-consistent analysis presented in section 4.1. However, if either /?o or 
A are non-zero then j3 is renormalized and the effect on n e is encoded in eqn. ( [4.2.14 ). 
In our case we set a R = cto = 1, Po = and A / 0. The equations ( 4.2.14|) , ([4.2.15 ) 
then give 



A 
2 







9k p 



(4.2.16) 



For small A and kq = we deduce that 



(— 

V187T 



1/6 



A l/3 



+ 



1 

18tt 



1/3 



A 2 / 3 . 



(4.2.17) 



The data for K e versus A for kq = is shown in figure and we see that for small 
A the simulation results are not compatible with A 2 / 3 behaviour. We shall see below 
that this is not rectified in the two-loop self-consistent calculation. However, in this 
one-loop calculation there is a considerable enhancement in the dependence of n e on 



A, whereas in the self-consistent calculation of section 4.1, in which the generation of 



the new vertex coupled to the vorticity u was not included, there is no effect at all at 
one-loop order and only a mild effect at two-loop order. The equations (4.2.16) can 
be solved numerically. For example, for kq = 0, A = 1 we find (3 R = 0.3456 and the 
effective velocity field is predicted to be 



u + (5 R u 



(4.2.18) 
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which clearly leads to an enhanced effective diffusivity, K e = 0.5207, compared with 
K e = 0.4090 from the two-loop calculation of the previous section. We believe that 
we have qualitatively captured the mechanism responsible for the enhancement of the 
effective diffusivity by helicity. 

The one-loop calculation is limited because it is not accurate at A = unlike the 
two-loop calculation. We have investigated the two-loop extension of the self-consistent 
approach when the new vertex with coupling (3 is included. This is more involved and 
the integrals were done numerically. We present the final results below. 

The two-loop self-consistent equations are 



where Bi and C\ are given in eqn. ( 4.2. 15| ) and B 2 and C 2 are evaluated numerically 
to be (or = Qo = 1) 

B 2 = -(0.0047A + 0.0095/3 + 0.0180/3A 2 + 0.0644/? 2 A + 0.0423/? 3 + 

0.0252/? 3 A 2 + 0.0287/? 4 A) , 
C 2 = 0.0088 - 0.0034A 2 - 0.0165/3A - 0.0315/3 2 - 0.0110/? 2 A 2 - 

0.0514/3 3 A + 0.0090/3 4 - 0.0337/3 4 A 2 . (4.2.20) 

Eqns. ( |4.2.19j ) can be rearranged to give 



Ak+^ + ±(C 2 -C 2 1 -B 1 ^) = 



-A/3+^ + ^(B 2 -2B 1 C 1 -B 1 ^) = 0. (4.2.21) 

K 2 e k% d[3 

We shall set [3q = from now on. These equations contain the accurate two-loop self- 
consistent result which fits the data for all kq at A = 0. For A / these equations 
are solved numerically and the results are compared with simulation data in figure 



- 11 As in the one- loop case the behaviour for small A is clearly incorrect and there 
is no quantitative agreement with the data. However, there is a clear enhancement 
of K e due to the inclusion of the vorticity vertex and associated coupling (5 and for 
A = 1.0 we find n e = 0.5959. This mechanism for enhancing turbulent diffusion cannot 
obviously be deduced from perturbative considerations. It arises from a resummation 
of diagrams which give an expression that analytically continues between the regions 
where A 2 <C Uq/kq/c 4 and A 2 S> "Uq/kq & 4 . The effective diffusivity in the former region 
is well predicted by perturbation theory, but not so for parameters in the latter region. 
Although quantitative agreement is not good the important point is that a mechanism 
has been discovered which gives a strong enhancement to the value of n e for non-zero A 
even in the one-loop approximation, whereas in the self-consistent theory for K e alone, 
discussed in the previous section, there is no effect of helicity at all on the value of 
K e at one-loop order. The obvious reason for the discrepancy in this new approach 
is that the approximations made are much too crude. A more refined calculation 
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would use a functional self-consistent method for Uj(k',k) (eqn. ( 4.2.3|) ) and Q(k 2 ) 
(eqn. (|4.2.4| )). Although a computationally formidable task, this is likely to encode the 
correct behaviour much more accurately than does our low-wavenumber approximation. 

The origin of the dip in figure || in the curve K e versus kq for A = 1 is unexplained 
by the theory presented so far. 



The Functional Hartree-Fock Method 



This approach goes some way towards including effects omitted in the low-wavenumber 
approximation. The version presented here is deficient in that the prediction for n e 
when kq = A = is not as accurate as the two-loop self-consistent approach but the 
advantage is that Q.(k 2 ), eqn. ( 4.2.4 ), is treated as function to be determined self- 
consistently by the Hartree-Fock equations. The vertices are still treated in the low- 
wavenumber approximation and, as in the previous section, they are parameterized by 
a and (3. 

The integral equation to be satisfied by Q(k 2 ) and the one-loop equation satisfied 
by the vertex function, which is the same as the one-loop self-consistent equation, are 
shown in figure ^. Note that, unlike the self-consistent calculation of the previous 
section, only one of the vertices in the one-loop self-energy is replaced with the full 
vertex since this gives the correct counting of diagrams when the equations are iterated. 
The self-consistent case is different because the augmented vertex is already present 
in the perturbation theory and corrections are implemented by counter-terms. The 
approximation for the vertices in eqn. ( 4.2.5| ) is used and j3 is determined using eqn. 
(|4.2.14|) : 

Bi 



= Po 



+ — 



(5.1) 



where B\ is given in eqn. ( |4.2.15| ) and using n e = kq + 17(0) . The Hartree-Fock 
equation to be satisfied by f2(fc 2 ) is then (after some reduction and using cj) in eqn. 

GED) 



Q(k 2 ) 



1 

F 



-k 2 /2 



3V27T 



dp 



(pk cosh(pk) — sinh(pfc)) e 



0\ 



d 3 p 
(2^p 



pk (K + n( P 2 )) 

|p + k| ^(|p + k|)(A;V-(k-p) 2 ) 



+ 



P 2 (k + Q( P 2 )) 



(5.2) 



These equations are solved by discretizing the wavevector, performing the integrals nu- 
merically and iterating the equations to converge to a solution for Q(k 2 ). The effective 
diffusivity is then K e = kq + 17(0) . 

The results are shown in figures ^ - 11, where it is clear that while the numerical 
value predicted at Ko = A = is not accurate, the behaviour for small A is more in 
keeping with the simulation results. This strengthens our belief that an analysis which 
treats the propagator and vertex as functions to be determined self-consistently is likely 
to reproduce the desired properties more accurately. However, at kq = 0, A = 1 the 
value predicted is K e = 0.5070, still much less than the simulation value of 0.8018(7). 

In principle, the vertex may be treated as a function in the same manner as Q(k 2 ). 
This is prohibitively expensive in memory and computer time but might be possible if 
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some simplification of the function form were implemented. Also, while the equation 
for Q(k 2 ) is already exact at one-loop, that for the vertex is not and we cannot preclude 
that higher loop corrections might be important. We have not pursued this approach. 

We note that in this approach, as with those of the previous sections, the marked 
dip 111 Hvg ELS cl function of kq for the larger values of A is not reproduced. 

6 The Renormalization Group 

In the previous section we presented an analysis based on the assumption that the 
large-scale advection is controlled by an effective transport equation dominated by 
the terms containing only one and two derivatives. This method is related to the 
renormalization group (RG) methods which have proved very successful in predicting 
exponents in critical phenomenon. In the RG approach a large wavenumber cutoff, 
A, is introduced and the advection on scales larger than L = 2tt/A is assumed to 
be described by an effective transport equation, in principle containing terms with an 
arbitrarily high number of derivatives. The parameters in this equation are functions 
of A in order to account for the effect of advection at the scales smaller than L which 
have been excised. In the limit A — ► the effective equation, by dimensional analysis, 
takes a simple form dominated by terms with few derivatives and with associated 
effective or "renormalized" parameters. In this way the effective equation takes a form 
similar to that used in the previous section. There is a difference, though, because 
any practical application of these schemes requires a drastic truncation of the operator 
space: especially in the RG method where it is impossible to compute the flow with 
changing A for very many parameters in the effective transport equation. Unlike the 
situation in critical phenomena there are no infra-red divergences in the theory and 
the notion of a relevant operator is not applicable. It is then a matter of trial and 
error to determine whether the approach used captures the vital features controlling 
the flow. The simplest renormalization scheme is to calculate the renormalization to 
the diffusivity k(A) and to the vertex associated with the coupling of the random field 
or externally applied drift. In the case of gradient flows we demonstrated in reference 
that this scheme yields exact results in one and two dimensions and an extremely 
accurate, although not exact, result in three dimensions. It is, in general, much harder 
to calculated the renormalized parameters such as n e than the associated exponents, 
and so success in |?J suggests that some insight may be gained using RG methods in 
other similar problems. 

In this section we present a RG calculation of K e . The vertex renormalization is done 
but multiple vertex renormalization is neglected which means that the renormalized 
velocity field remains Gaussian. Consequently, after integrating out the random field 
down to wave number A we postulate that the equation for the effective Green function 
can be approximated, for all A, by an equation of the same form as the original one 
(eqn. P): 

k(A)V 2 G(x,A) - u(x,A) • VG(x,A) = -<5(x) , (6.1) 

where /c(A) is the running renormalized diffusion constant and ua is the renormalized 
velocity field. Since we renormalize the vertex functionally the field correlation function 
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will flow under the RG as 



ft(k,AK(K'.A» - { «+ k '>^- A > W < A (62) 

where 

Eij(k,A) = &(k,A)(k 2 5ij - kikj) + tf(k, A)ie imj k m . (6.3) 

One finds that the renormalized field is still incompressible. We shall compute the flow 
equations for «(A), <£(k, A) and ^(k, A) as A varies. 

The change in «(A) on integrating out wave vectors in the shell (A, A — 5 A) is 

*«(A) = -^^A 2 <S>(AA)5\ (6.4) 

If one calculates the vertex renormalization and treats it as an addition to the 
random field one finds, 

(2vr) 3 7 A _5A K 2 (A)g' 2 (q + q') 2 

where we have only kept the vertex term to lowest order in k as it is only this term that 
contributes to the one loop diffusivity renormalization. In addition if one assumes that 
it is the low wave number (long distance) renormalization of the velocity field which is 
important for the effective diffusivity one finds to lowest order in q, 

$Ui(q, A) « - 1 qj u k ( q ,A) [ A Fij ^ Wk H (6-6) 
(27r) d ^(A) J as A <? 4 

In the case where no helicity is present we see that the vertex is not renormalized. 
However (as pointed out in the section on perturbation theory) whilst the helicity does 
not contribute to the diffusivity renormalization at a one loop level it renormalizes the 
velocity field. The renormalization is zero at order in q but has an order 1 effect: 

5ui(q,A) » - 67r 2 K 2( A ) tijkQjUiiq, A)¥(A, A)8A (6.7) 

In real space therefore, the renormalization is of the form 

u -> u + AA/3(A)V x u. (6.8) 

Using the renormalization of u one may compute the flow of and thus and \P to 
obtain the one loop functional RG equations: 

° H ' -A 2 <D(A,A) 



dA 3tt 2 k(A)' 
8A 3tt 2 « 2 (A) W ' ;n ' j 

The integration of the eqns. (|6.9| ) is from A = oo to with the initial conditions 

/c(oo) = Kq 
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*(g,oo) = 

*(g,oo) = *(g). (6.10) 

When there is no helicity there is no vertex renormalization at the order we are 
considering and therefore we may integrate the equations directly to obtain 

Ke = (/«g+2u§/9)3 . (6.11) 

This agrees with the one loop perturbation result, as it should, and in the case kq = 
we find that K e = \/2/3 = 0.47140 which is quantitatively not very close to the 
numerically measured result, K e = 0.3697. However, the discrepancy is sensitive to the 
form assumed for the effective diffusion equation. In our case this is given by eqn. ( |6.lD 
which is clearly inadequate since u is not renormalized when A = 0. An improvement 
can only be made by including terms with higher derivatives of u. This is similar to 
parameterizing the non-helical form factor V R of eqn. ( [4.2.3 ) with a function of external 



momenta rather than approximating it by a constant, or, which is not renormalized. 
This is a possible avenue of research but we have not yet followed it. 

In contrast, for A 7^ 0, u is renormalized and the effect on K e is significant because the 
helical form factor W R , eqn. ( 4.2.3Q is renormalized at low wavenumber as parametrized 



by /3(A) above. The RG equations may be integrated numerically and is compared with 
simulation in figures |j] - 11. Although the results are not quantitatively accurate, they 



capture the qualitative behaviour seen in the simulations. In particular, the RG predicts 
the large enhancement as a function of A seen in the data and also predicts the dip 
observed in the graph of K e versus kq for sufficiently large A. 

Indeed, the qualitative success of the method suggests that the difficulty in obtaining 
predictions that are more accurate might lie with the inadequacy of the simple ansatz 
when applied to the case when A = 0. The effect of helicity is nevertheless well captured 
in this approach because this effect is dominated by the renormalization of /3(A). 

A technical point in the numerical integration is that kq = =^ k(oo) = 0, and the 
evolution equations are ill-defined in the limit A — > 00. This problem is easily rectified 
by making kq very small but non-zero. The integration procedure is then well-defined 
and the results are insensitive to the exact value of kq in this case. 

We therefore believe that although the renormalization procedure is not quanti- 
tatively accurate (as should be expected as it does not give very accurate results in 
the absence of helicity), it successfully incorporates the underlying mechanism for the 
enhancement of the diffusivity by helicity at small bare molecular diffusivity. 

7 Discussion and Conclusions 

In this paper we have studied the problem of turbulent advection of a scalar field by an 
incompressible flow with helicity A, < A < 1.0, and background molecular diffusivity, 
kq. We have performed computer simulations of the advection for flows with properties 



described in eqns. ( |2,lD to (|2.4j ), and compared the long-range effective parameters 
describing the time evolution of the scalar field with various schemes of calculation. In 
particular, we have concentrated on how the effective diffusivity, K e , depends on kq and 
A. In earlier work we found an strong anomalous enhancement of function of 

A for kq = 0.0 0] which was unexplained theoretically, and this is the motivation for 
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the present study. In that earlier work the turbulent velocity field was time dependent 
whereas here it is not. This allows for easier calculation whilst still reproducing the 
effect. 

The important region for discussion can be seen from the simulation data, figures 
- [□], to be k < 0.2; for lar ger Ato the molecular diffusivity begins to dominate and 
not only is the effect of helicity suppressed but also the many schemes of calculation 
give good approximations for K e . For A = 0.0 we find that the two-loop self-consistent 
calculation of n e reproduces the simulation data for all kq very closely indeed, as is 
seen in figure |7| and described in section ||[ The other schemes also plotted are much 
less accurate in the region of interest. Ordinary perturbation theory is not convergent 
in this region and will be discussed no further. The reason why the Hartree-Fock and 
(RG) methods are less accurate is that the vertex function V(k',k), eqn. (4.2.3), is not 
renormalized for low wavenumber which means that the associated coupling a is not 
renormalized. The Hartree-Fock method at A = sums rainbow diagrams but does not 
include any diagrams corresponding to a vertex correction, unlike the self-consistent 
theory. In the self-consistent theory the one-loop prediction for kq = is n e = 1/3 and 
the two-loop terms modify this by 5n e ~ 0.04, of which the two-loop cross diagram 
in figure |3| contributes only 10%, or SK e ~ 0.004. In omitting similar terms to this 
latter one the Hartree-Fock approximation should therefore not be expected to be too 
discrepant and this is seen to be the case. The RG calculation gives a form which 
must yield the simple one-loop perturbation theory expression at large kq but allows 
continuation to kq = 0; this is given in eqn. ( |6.11 ). In the case of gradient flows the 
RG approach is remarkably successful Q and this is attributed to the fact that in that 
case the primitive vertex is renormalized at low wavenumber. 

The reason for examining schemes alternative to self-consistent methods is that 
for A > agreement between simulation data and theory is poor and it is necessary 
to investigate different approaches in order to test different hypotheses for a simple 
description of the observed anomalous effect. 

For kq > 0.2 all schemes except ordinary perturbation theory begin to show rea- 
sonable agreement with the data, and for kq > 0.5 all schemes clearly reproduce the 
results. We concentrate on results for kq < 0.2 and the anomalous enhancement of n e 
by helicity in this region is seen in figure |K] where K e is plotted against A for kq = 0, 
and is characterized by a rapid rise in K e for A > 0.2 . An alternative aspect is seen 
in figure ^ where K e is plotted against kq for A = 1.0. The significant dip at kq = 0.2 
is due to the large effect of helicity on n e smaller kq compared with the much reduced 
effect at kq ~ 0.2. It has proved very difficult to convincingly explain these features. 
However, we have been able to suggest mechanisms which show the presence of helicity 
in the flow can produce a large change in K e from the non-helical value, and even though 
these have not yet yielded quantitative predictions they do point towards a reasonable 
explanation. 

The basic idea is to recognize that the effective equation governing the advection 
should contain terms not present in the original equation. The terms can be thought 
of as being induced in the low- wavenumber effective theory by integrating out higher 
wavenumbers. This may also be viewed as the renormalization of the related vertex 
functions of the theory corresponding to a selective resummation of diagrams. In our 
approach we have assumed that a low- wavenumber approximation will be valid and so 
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such terms will contain the minimum number of derivatives. These ideas can be imple- 
mented in different ways and we tried self-consistent, Hartree-Fock and renormalization 
group approaches. The self-consistent and Hartree-Fock methods are based on the ef- 



fective evolution equation ( 4.2.1 ) which corresponds to a low-wavenumber enhancement 
of the flow velocity ur = u + /3rw, where uj is the vorticity. We performed a two-loop 
calculation self-consistent in both K e and (3r as described in eqn. (14.2.19[) and figure |. 



The results show that a strong enhancement in K e is predicted but that the magnitude 
for A = 1.0, kq = 0.0 is too small and the form of the dependence of n e on A disagrees 
with the data. This is particularly true at small A for kq = 0.0 where from figure [H] we 
see that K e is only weakly dependent on A, A < 0.3 whereas we predict K e ~ a + bX p 
for fractional p: the one-loop result is p = 2/3. 

The Hartree-Fock method, shown diagramatically in figure [6|, computes the com- 
plete propagator in terms of Q(k 2 ) (eqn. (|4.2.4 )) as a sum of the rainbow diagrams 



generated from the effective equation ( f4.2.1| ) with fin given by the one-loop result 
(3r = —B\{(5r) / K 2 where B\{[3r) is given in eqn. ( 4.2.15 ) and K e = kq + $7(0). A rise 



in n e with A is predicted but not one big enough to agree with the data. However, the 
behaviour at small A is a much slower rise which is more in keeping with the data than 
the self-consistent prediction in this region. 

The renormalization group method is a different approach in that it considers a 
running diffusivity k(A) and velocity field u(x, A) which satisfy k(oo) = kq, k(0) = 
K e , u(x, A) = u(x). There are three coupled RG flow equations, fl6,9|), for k(A) and the 
two running spectral functions <£((/, A), *$>(q, A) which correspond to the definitions in 
eqn. ( |6 . 3| ) . The important feature of the RG flow equations is that k(A) appears in the 
denominators. The numerators are suppressed at large A by the spectral functions and 
so the major contribution is from intermediate values of A: A ~ k$. This contribution 
is strongly enhanced for small kq and results in the prediction of the dip structure 
observed in the data but not predicted by the other methods. 

^From our investigation we must expect that a proper explanation of the observed 
effects will require the correct effective equation and the consequent generation of new 
vertices but that unlike the A = case the low-wavenumber approximation will be 
insufficient since although an enhancement is predicted for A ^ the rapid rise is not 
reproduced and no dip is observed. The RG method suggests that the main contribution 
is from wavenumbers k ~ k$ supporting this latter conclusion. A successful approach 
should therefore include more terms in the effective flow equation in combination with 
an RG approach. The challenge is to obtain accurate results for all A including A = 
by such a technique. Work in this direction is currently underway. 
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Figure 1: The vertices occurring in the perturbation schemes: (a) the primitive velocity held 
vertex; (b) the primitive vorticity vertex; (c) the bare complete vertex of the effective diffusion 
equation; (d) the renormalized complete vertex approximated as a sum of renormalized vertices 
associated with the velocity held and the vorticity. 
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Figure 2: The graphs that contribute to two-loop simple perturbation theory for G(k). 
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Figure 3: The perturbation expansion to two-loops of the self-consistent relation for K e 
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Figure 4: An example of the kind of vertex graph that must be evaluated in the solution of 
equations shown in figure]^ once the approximation for \J R given in eqn. (4.2.5) and shown in 
figure |l| has been used. This graph is labelled Tp aa . 
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Figure 5: The perturbation expansion to two-loops of the general self-consistent condition re- 
lating the self-energy, £(£;), and vertex, XJ R (k, k'), functions. The vertex function is represented 
by the circle with inset 'R', full Green function, G(k), by the filled box and AU = XJ R — U°. 
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Figure 6: The integral Hartree-Fock equations for f2(k) and XJ R in terms of the general vertices 
U° and . The full Green function, G(k), is denoted by the filled box. The equation for 
f2(k) is correct to all loop orders but the equation for XJ R is correct to one-loop order only. 
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Figure 7: n e versus kO for fixed helicity A = 0.0. The simulation data are shown (o) to 
be compared with the predictions of two- loop self-consistent perturbation theory (solid), the 
Hartree-Fock calculation (long-dashed), the renormalization group (dashed) and ordinary per- 
turbation theory (dot-dashed) 
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Figure 8: n e versus k0 for fixed helicity A = 0.4. The simulation data are shown (o) to 
be compared with the predictions of two-loop self-consistent perturbation theory in n e (solid) 
and in « e ,/3 (dotted), the Hartree-Fock calculation (long-dashed), the renormalization group 
(dashed), and ordinary perturbation theory (dot-dashed) 
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Figure 9: n e versus k0 for fixed helicity A = 1.0. The simulation data are shown (o) to 
be compared with the predictions of two-loop self-consistent perturbation theory in n e (solid) 
and in « e ,/3 (dotted), the Hartree-Fock calculation (long-dashed), the renormalization group 
(dashed), and ordinary perturbation theory (dot-dashed) 
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Figure 10: K e versus A for fixed diffusivity kq — 0.0. The simulation data are shown (o) to 
be compared with the predictions of two-loop self-consistent perturbation theory in n e (solid) 
and in K e ,f3 (dotted), the Hartree-Fock calculation (long-dashed), the renormalization group 
(dashed) 
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Figure 11: « e versus A for fixed diffusivity kq = 0.2. The simulation data are shown (o) to 
be compared with the predictions of two-loop self-consistent perturbation theory in n e (solid) 
and in K e ,f3 (dotted), the Hartree-Fock calculation (long-dashed), the renormalization group 
(dashed) 
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